Fractional-Period Excitations in Continuum Periodic Systems 
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We investigate the generation of fractional-period states in continuum periodic systems. As an 
example, we consider a Bose-Einstein condensate confined in an optical-lattice potential. We show 
that when the potential is turned on non-adiabatically, the system explores a number of transient 
states whose periodicity is a fraction of that of the lattice. We illustrate the origin of fractional- 
period states analytically by treating them as resonant states of a parametrically forced Duffing 
oscillator and discuss their transient nature and potential observability. 
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I. INTRODUCTION 

In the past few years, there has been considerable inter- 
est in both genuinely discrete and continuum but periodic 
systems pj. These arise in diverse physical contexts 0, 
including coupled waveguide arrays and photorefractive 
crystals in nonlinear optics , Bose-Einstein condensates 
(BECs) in optical lattices (OLs) in atomic physics Q, 
and DNA double-strand dynamics in biophysics 5] . One 
of the most interesting themes that emerges in this con- 
text is the concept of "effective discreteness" induced by 
continuum periodic dynamics. There have been many 
efforts both to derive discrete systems that emulate the 
dynamics of continuum periodic ones and to obtain 
continuum systems that mimic properties of discrete ones 
Q. Additionally, the connection between discrete and 
continuum systems in various settings is one of the main 
research thrusts that has emerged from studies of the 
Fermi-Pasta- Ulam problem • 

This paper examines a type of excitation, not previ- 
ously analyzed (to the best of our knowledge), with the 
intriguing characteristic that it can be observed in con- 
tinuum periodic systems but cannot be captured using a 
genuinely discrete description of the same problem. The 
reason for this is that these states bear the unusual fea- 
ture that their length scale is a fraction of that of the con- 
tinuum periodic potential. Thus, the "fractional-period" 
states reported in this paper are not stationary states of 
the latter problem, but rather transient excitations that 
persist for finite, observable times. 

To illustrate these fractional states, we consider the 
example of a trapped BEC in which an OL potential is 
turned on (as a non-adiabatic perturbation) |4( . Our re- 
sults can also be applied in the context of optics by con- 
sidering, for example, the effect of abruptly turning on 
an ordinary polarization beam in a photorefractive crys- 
tal Q. Our particular interest in BECs is motivated by 
recent experiments jg, where after loading the conden- 
sate in an OL, the amplitude of the pertinent standing 



wave was modulated and the resulting excitations were 
probed. These findings were subsequently analyzed in 
the framework of the Gross-Pitaevskii equation in 0, 
where it was argued that a parametric resonance occurs 
due to the OL amplitude modulation. These results were 
further enforced by the analysis of ^3] ■> which included a 
computation of the relevant stability diagram and growth 
rates of parametrically unstable modes. The results of @ 
were also examined in by treating the Bose g as as a 
Tomonaga-Luttinger liquid. 

A similar experiment, illustrating the controllability of 
such OLs, was recently reported in , where instead of 
modulating the amplitude of the lattice, its location was 
translated (shaken) periodically. This resulted in mixing 
between vibrational levels and the observation of period- 
doubled states. Such states were predicted earlier p^3.ll4| 
in both lattice (discrete nonlinear Schrodinger) and con- 
tinuum (Gross-Pitaevskii) frameworks in connection to 
a modulational instability 0, ^| and were also exam- 
ined recently in |17| . Period-multiplied states may exist 
as stationary (often unstable) solutions of such nonlin- 
ear problems and can usually be captured in the relevant 
lattice models. 

To obtain fractional-period states, which cannot be 
constructed using Bloch's theorem [lj|, we will consider 
a setting similar to that of 0, akin to the experiments 
of |8(- However, contrary to the aforementioned earlier 
works (but still within the realm of the experimentally 
available possibilities of, e.g., Ref. @), we propose ap- 
plying a strong non-adiabatic perturbation to the sys- 
tem (which originally consists of a magnetically confined 
BEC) by abruptly switching on an OL potential. As 
a result, the BEC is far from its desired ground state. 
Because of these nonequilibrium conditions, the system 
"wanders" in configuration space while trying to achieve 
its energetically desired state. In this process, we monitor 
the fractional-period states as observable transient exci- 
tations and report their signature in Fourier space. After 
presenting the relevant setup, we give an analysis of half- 
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period and quarter-period states in a simplified setting. 
We illustrate how these states emerge, respectively, as 
harmonic and 1:2 super harmonic resonances of a para- 
metrically forced Duffing oscillator describing the spatial 
dynamics of BEC standing waves (see Appendices A and 
B for details). We subsequently monitor these states in 
appropriately crafted numerical experiments and exam- 
ine their dependence on system parameters. Finally, we 
also suggest possible means for observing the relevant 
states experimentally. 

The rest of this paper is organized as follows. In Sec- 
tion II, we present the model and the analytical results. 
(The details of the derivation of these results are pre- 
sented in appendices; Appendix A discusses half-period 
states and Appendix B discusses quarter-period states.) 
We present our numerical results in Section III and sum- 
marize our findings and present our conclusions in Sec- 
tion IV. 



II. MODEL AND ANALYSIS 
A. Setup 

A quasi-lD BEC is described by the dimensionless 
Gross-Pitaevskii (GP) equation pUli^. 

i ^=-\^+9m+v( X ,m, a) 

where t/>(x, t) is the mean-field wavefunction (with atomic 
density \ip\ 2 rescaled by the peak density n ), x is mea- 
sured in units of the healing length £ = h/ ^noguym 
(where m is the atomic mass), t is measured in units 
of £/c (where c = \fnoguj/rn is the Bogoliubov speed 
of sound), g\D = 2hio±a is the effective ID interaction 
strength, u>± is the transverse confinement frequency, a 
is the scattering length, and energy is measured in units 
of the chemical potential /i = giDn-o- The nonlinearity 
strength g (proportional to a) is taken to be positive in 
connection to the 87 Rb experiments of [j| . The potential, 

V(x,t) = i fiV + V H (t)[l + A sin (ujt)] sin 2 (qx), (2) 

consists of a harmonic (magnetic) trap of strength O = 
fv^>xl 'giD^o (where u> x is the longitudinal confinement fre- 
quency) and an OL of wavenumber q, which is turned on 
abruptly at t — [via the Heaviside function H(t)}. The 
lattice depth, given by Vo [1 4- Asm (cut)], is periodically 
modulated with frequency lo. 

Before the OL is turned on (i.e., for t < 0), the mag- 
netically trapped condensate is equilibrated in its ground 
state, which can be approximated reasonably well by the 
Thomas-Fermi (TF) cloud utf = ^/max{0, /Jo — V(x)}, 
where /io is the normalized chemical potential [3j. The 
OL is then abruptly turned on and can be modulated 
weakly or strongly (by varying A) and slowly or rapidly 
(by varying u). 



To estimate the physical values of the parameters in- 
volved in this setting, we assume (for fixed values of the 
trap strength and normalized chemical potential, given 
by f2 = 0.01 and /iq — 1, respectively) a magnetic trap 
with ujj_ = 2vr x 1000 Hz. Then, for a 87 Rb ( 23 Na) con- 
densate with ID peak density 5 x 10 7 m _1 and longitudi- 
nal confinement frequency lo x = 2ir x 6 Hz (2tt x 2.8 Hz), 
the space and time units are 0.4/im (1.25^m) and 0.27 
ms (0.57 ms), respectively, and the number of atoms (for 
g= 1) is Nm 4200 (12000). 

B. Analytical Results 

To provide an analytical description of fractional- 
period states, we initially consider the case of a homoge- 
neous, untrapped condensate in a time-independent lat- 
tice (i.e., 51 = A = 0). We then apply a standing wave 
ansatz to Eq. (JIJ to obtain a parametrically forced Duff- 
ing oscillator (i.e., a cubic nonlinear Mathieu equation) 
describing the wavefunction's spatial dynamics. As ex- 
amples, we analyze both half-period and quarter-period 
states. We discuss their construction briefly in the 
present section and provide further details in Appendices 
A and B, respectively. 

We insert the standing wave ansatz 

ip(x, t) = R(x) exp (-ifiot) exp [-i(V /2)t] (3) 

into Eq. to obtain 

R" + 5R + &R 3 + eV Rcos(Kx) = , (4) 

where primes denote differentiation with respect to x, 
S = 2//q, a = —2g, eVq = Vq, and k = 2q. 

We construct fractional-period states using a multiple- 
scale perturbation expansion pol |. defining rj = ex and 
£ = bx = (1 +ebi + e 2 &2 + • • • )x for stretching parameters 
bj. We then expand the wavefunction amplitude R in a 
power series, 

R = R +eR 1 + e 2 R 2 + e 3 R 3 + 0(e 4 ) . (5) 

Note that although £ and r\ both depend on the vari- 
able x, the prefactor e in 7/ indicates that it varies much 
more slowly than £ so that the two variables describe phe- 
nomena on different spatial scales. In proceeding with a 
perturbative analysis, we treat £ and rj as if they were 
independent variables (as discussed in detail in Ref. pojl 
in order to isolate the dynamics arising at different scales 
i2(J . We also incorporate a detuning into the proce- 
dure (in anticipation of our construction of resonant so- 
lutions) by also stretching the spatial dependence in the 
OL, which gives W(£) = eVo cos(k£) for the last term in 
Eq. |@J We insert Eq. J5J into Eq. (gj, expand the 
resulting ordinary differential equation (ODE) in a power 
series in e, and equate the coefficents of like powers of e. 

At each 0(s^), this yields a linear ODE in £ that Rj = 
Rj(£,rf) must satisfy: 

L^Rj] ee *y?L + SRj = hj(£, R k ,D l R k ) , (6) 
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where hj depends explicitly on £ and on Rk and its 
derivatives (with respect to both £ and 77) for all k < j. 
We use the notation D l Rk in the right hand side of 
Eq. 10 to indicate its functional dependence on deriva- 
tives of i?fe. In particular, because of the second deriva- 
tive term in Eq. these terms are of the form , 
, and . [See, for example, Eq. l(TT|) in Appendix 

A.] 

We scale a (see the discussion below) to include at least 
one power of e in the nonlinearity coefficient in order 
to obtain an unforced harmonic oscillator when j = 
(so that ho vanishes identically) 21]. At each order, we 
expand hj in terms of its constituent harmonics, equate 
the coefficients of the independent secular terms to zero, 
and solve the resulting equations to obtain expressions 
for each of the Rj in turn. (The forcing terms hj and 
the solutions Rj are given in Appendix A for half-period 
states and Appendix B for quarter-period states.) Each 
Rj depends on the variable 77 through the integration 
constants obtained by integrating Eq. |JB} with respect to 
£. The result of this analysis is an initial wavefunction, 
ip(x,0) = R(x), given by Eq. (J5J. 

We obtain half-period states of Eq. @ (and hence of 
the GP equation) by constructing solutions in harmonic 
(1:1) resonance with the OL (i.e., \/6 — k — 2q) [23| . 
To perform the (second-order) multiple-scale analysis for 
this construction (see Appendix A), it is necessary to 
scale the nonlinearity to be of size 0(e) (i.e., a — sa), 
where e is a formal small parameter. [The OL is also of 
size 0(e).] We show below that full numerical simulations 
of the GP equation with a stationary OL using initial 
conditions obtained from the multiple-scale analysis yield 
stable half-period solutions even for large nonlincarities. 
The oscillations in time about this state are just larger 
because of the O(l) nonlinearity. 

We also obtain quarter-period states of Eq. J3J (and 
hence of the GP equation) by constructing solutions in 
1:2 superharmonic resonance with the OL (i.e., y5 = 
2k = 4q). Because the 1:2 superharmonic resonant so- 
lutions of the linearization of (0J [that is, of the linear 
Mathieu equation] are 4th-order Mathieu functions |24| . 
we must use a fourth-order multiple-scale expansion (see 
Appendix B) to obtain such solutions in the nonlinear 
problem when starting from trigonometric functions at 
0(1). Accordingly, it is necessary to scale the nonlin- 
earity to be of size 0(e 4 ) (i.e., a = e 4 a). [The OL 
is still of size 0(e).] Nevertheless, as with half-period 
states, we show below that full numerical simulations of 
the GP equation with a stationary OL using initial con- 
ditions obtained from the multiple-scale analysis yield 
stable quarter-period solutions even for large nonlinear- 
ities. The oscillations in time about this state are again 
larger because of the 0(1) nonlinearity. 
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FIG. 1: (Color online) Evolution of the half-period state found 
by numerically integrating Eq. l|T). (Left) Space-time plot. 
(Right) Snapshots of the density \ip\ 2 at t = 0, t = 100, 
t = 200, and t = 300. The parameter values are f2 = A = 
(i.e., a time-independent OL and no magnetic trap), g = 1, 
V = l,q = n/8 (so VS = tt/4), e = 0.05, and 61 = 62 = -1. 
The dashed curve shows the OL potential. 



III. NUMERICAL RESULTS 

Having shown the origin of fractional-period states an- 
alytically, we now use numerical simulations to illustrate 
their dynamical relevance. 

First, we consider the case with a stationary OL (A = 
0) in the absence of the magnetic trap (Q = 0). We ex- 
amine the time-evolution of half-period states by numer- 
ically integrating Eq. JQ| with the initial condition given 
by Eqs. and (J5J. We show an example in Fig.^ where 
g = 1, Vq = 1, q = 7r/8, and e = 0.05. The half-period 
state persists for long times (beyond t — 300). The pa- 
rameter values correspond to a 87 Rb ( 23 Na) BEC with a 
ID peak density of 5 x 10 7 m" 1 confined in a trap with 
frequencies lo±_ = 2ix x 1000 Hz, lu x = 0, and number 
of atoms N « 4000 (12500). In real units, t = 300 is 
about 80 ms (170 ms). We similarly examine the time- 
evolution of quarter-period states by numerically inte- 
grating Eq. with the initial condition again given by 
Eqs. © and JSJ, but now for the case of the superhar- 
monic 1:2 resonance. In Fig. [21 we show an example for a 
similar choice of parameters as for half-period states [25| . 

We subsequently examine the generation of such states 
through direct numerical experiments using Eq. in 
the presence of magnetic trapping. Initially, we include 
the parabolic and periodic components of the potential 
but leave the potential time- independent, setting A = 0, 
/Jo = 1, Vq = 1, and q = v2/2. First, we consider 
the case of a weak parabolic trap with fl = 0.001. Us- 
ing the nonlinearity coefficient g = 1, we perform the 
numerical experiment as follows: We integrate the GP 
equation in imaginary time to find the "exact" ground 
state (in the absence of the OL) and then we switch on 
the OL at t = 0. We monitor the density \ip(x, i)| 2 , its 
Fourier transform \ ^(k, t)\ 2 , and the spectral components 
at k = 2q (the OL wavenumber) and k — 4q (half the 
wavenumber). Note that the spectrum also contains a 
"DC" component (at k = 0, corresponding to the ground 
state) as well as (very weak) higher harmonics. 

In general, the spectral component at k = 2q is sig- 
nificantly stronger than that at k = 4g (see the bottom 
panel of Fig. , implying that the preferable scale (pe- 
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FIG. 2: (Color online) Evolution of the quarter-period state 
found by numerically integrating Eq. Q . (As in Fig. Q the 
lattice is time-independent and the magnetic trap is not in- 
cluded.) We show snapshots of the density \ip\ 2 at t = 0, 



t = 100, t = 200, and t = 
Sl = A = 0,g = l,Vo = l,q 
and bi = 62 = &3 = 64 = — 1. 
potential. 



300. The parameter values are 
= tt/8 (so VS = tt/2), e = 0.05, 
The dashed curve shows the OL 
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FIG. 3: (Color online) Top panels: Snapshots of IVK^jt)! 2 
(left) and its Fourier transform |$(fc,f)| 2 (right) for the case 
of a time-independent lattice and a magnetic trap at times t = 
1.2, t = 7.3, and t = 7.7. The parameter values are f2 = 0.001, 
^0 = 1, Vo = 1, q = y/2/2, and A = 0. The dashed curves 
(in the left panels) show the OL. Bottom panel: Evolution of 
|*(4g)| 2 (thin curve) and |#(2<?)| 2 (bold curve). For t = 1.2, 
the density has the same period as the OL (i.e., k — 2q). 
For t = 7.3 (t — 7.7), we observe the formation of a quasi- 
harmonic (non-harmonic) half-period state with wavenumber 
k = 4q. 



riod) of the system is set by the OL. This behavior is most 
prominent at certain times (e.g., at t = 1.2), where the 
spectral component at k — 2q is much stronger than the 
other harmonics. Nevertheless, there are specific time in- 
tervals (of length denoted by r) with |*(4(j)| 2 > |<F(2g)| 2 , 
where we observe the formation of what we will hence- 
forth call a "quasi-harmonic" half-period state. For ex- 
ample, one can see such a state at t = 7.3. The purpose of 
the term "quasi-harmonic" is to characterize half-period 



states whose second harmonic (at k = Aq, in this case) 
is stronger than their first harmonic (at k = 2q). As 
mentioned above, such states have an almost sinusoidal 
shape, like the wavefunctions we constructed analytically. 

One can use the time-evolution of the spectral com- 
ponents as a quantitative method to identify the for- 
mation of half-period states. This diagnostic tool also 
reveals a "revival" of the state, which disappears and 
then reappears a number of times before vanishing com- 
pletely. Furthermore, we observe that other states that 
can also be characterized as half-period ones (which tend 
to have longer lifetimes than quasi-harmonic states) are 
also formed during the time-evolution, as shown in Fig. 
3 at t = 7.7. These states, which we will hereafter call 
"non-harmonic" half-period states, have a shape which 
is definititvely non-sinusoidal (in contrast to the quasi- 
harmonic states); they are nevertheless periodic struc- 
tures of period k = Aq. In fact, the primary Fourier peak 
of the non-harmonic half-period states is always greater 
than the secondary one. Such states can be observed for 
times t such that the empirically selected condition of 
|*(2g,i)| 2 < 3\^{Aq,t)\ 2 is satisfied. 

We next consider a stronger parabolic trap, setting 
fl = 0.01. Because the system is generally less homo- 
geneous in this case, we expect that the analytical pre- 
diction (valid for Q = 0) may no longer be valid and that 
half-period states may cease to exist. We confirmed this 
numerically for the quasi-harmonic half-period states. 
However, non-harmonic half-period states do still appear. 
The time-evolution of the spectral components at k = 2q 
and k = Aq is much more complicated and less efficient as 
a diagnostic tool, as \^{Aq,t)\ 2 < \^(2q,t)\ 2 for all t. In- 
terestingly, the non-harmonic half-period states seem to 
persist as is increased, even when the resonance condi- 
tion \/S = k = 2q is violated. For example, we found that 
for time-independent lattices (i.e., A = 0), the lifetime of 
a half-period state in the resonant case with q = y/2/2 
(recall that /io = 1) is r w 1.72, whereas for q = 1/2 it is 
r w 0.84. Moreover, the simulations show that the life- 
times become longer for periodically modulated OLs (us- 
ing, e.g., A=l; also see the discussion below). In partic- 
ular, in the aforementioned resonant (non-resonant) case 
with q = v2/2 (q = 1/2), the lifetime of the half-period 
states has a maximum value, at lu — 1.59 (u> — 0.75), 
of r w 8.24, or 4.7 ms (r » 5.72, or 3.3 ms) for a 23 Na 
condensate. We show the formation of these states in the 
top panels of Fig. ^ 

We also considered other fractional states. For exam- 
ple, using the same parameter values as before except for 
yS = 2k (so that q = y/2/A), we observed quarter-period 
transient states with lifetime r w 2.9. These states oc- 
curred even in the non-resonant case with q = 1/4 (yield- 
ing 7 « 1). We show these cases (for a temporally mod- 
ulated lattice with modulation amplitude A = 1) in the 
bottom panels of Fig. 0] In Fig. |SJ we show the lifetime 
r for the half-period and quarter-period states as a func- 
tion of A. Observe that the lifetime becomes maximal 
(for values of A < 1.5) around the value A = 1 considered 
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FIG. 6: (Color online) Density profiles for a fractional state 
with y/S = (5/4)k, so that q = 2\/2/5 (for /io = 1, lattice- 
modulation amplitude A = 1, and lattice-modulation fre- 
quency w = 1.59). The system oscillates in time between half- 
period (top panel) and third-period (bottom panel) states. 



FIG. 4: Half-period (top panels) and quarter-period (bottom 
panels) states for the case of a temporally modulated OL. 
The state in the top left panel has wavenumber q = \/2/2 
(resonant), and that in the top right panel has wavenumber 
q — 1/2 (non-resonant). Similarly, the states in the bottom 
panels have wavenumbers q — %/2/4 (left, resonant) and q = 
1/4 (right, non-resonant). 
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FIG. 5: The lifetime r of the half-period (upper curve) 
and quarter-period (lower curve) states as a function of the 
lattice-modulation amplitude A. For the half-period (quarter- 
period) state, we show the case with wavenumber q = V2/2 
(q = \/2/4). In both examples, we used /io = 1 and lattice- 
modulation frequency to = 1.59. 



above. Understanding the shape of these curves and the 
optimal lifetime dependence on A in greater detail might 
be an interesting topic for further study. 

Finally, we also examined non-integer excitations, for 
which the system oscillates between the closest integer 
harmonics. For example, in Fig. we show a state cor- 
responding to \fl = (5/4)k, so that q = 2^/2/5 (with 
= !)■ The system oscillates between the k = Aq 
(half-period) and k = 6q (third-period) states. Recall 
that the case presented in the top right panel of Fig. 0] 
(with q = 1/2) was identified as a "non-resonant half- 
period state" (the resonant state satisfies q — v2/2). 
Here it is worth remarking that this value, q = 1/2, is 



closer to q = \/2/3 (characterizing the third-period state) 
than to s/2/2. Nevertheless, no matter which charac- 
terization one uses, the salient feature is that the value 
q = 1/2 is nonresonant and lies between the third-period 
and half-period wavenumbers. Accordingly, the respec- 
tive state oscillates between third-period and half-period 
states. Thus, in the case shown in the top right panel of 
Fig. 21 the third-period state also occurs (though we do 
not show it in the figure) and has a lifetime of r = 4.2 
(2.4 ms), while for q = 2^/2/5 its lifetime is r = 2.92 
(1.67 ms). This indicates that states with wavenumbers 
closer to the value of the resonant third-period state have 
larger lifetimes. This alternating oscillation between the 
nearest resonant period states is a typical feature that we 
have observed for the non-resonant cases. 



IV. CONCLUSIONS 

In summary, we have investigated the formation and 
time-evolution of fractional-period states in continuum 
periodic systems. Although our analysis was based on a 
Gross-Pitaevskii equation describing Bose-Einstein con- 
densates confined in optical lattices, it can also be ap- 
plied to several other systems (such as photonic crystals 
in nonlinear optics). We have shown analytically and 
demonstrated numerically the formation of fractional- 
period states and found that they may persist for suffi- 
ciently long times to be observed in experiments. The 
most natural signature of the presence of such states 
should be available by monitoring the Fourier transform 
of the wavepacket through the existence of appropriate 
harmonics corresponding to the fractional-period states 
(e.g., k = Aq for half-period states, k — 8q for quarter- 
period states, etc.). 

It would be interesting to expand the study of the 
parametric excitation of such states in order to better 
understand how to optimally select the relevant driving 
amplitude. Similarly, it would be valuable to examine 
more quantitatively the features of the ensuing states as 
a function of the frequency of the parametric drive and 
the parabolic potential. 
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V. APPENDIX A: ANALYTICAL 
CONSTRUCTION OF HALF-PERIOD STATES 

To construct half-period states, we use the resonance 
relation ^/S = k and the scaling a — ea, so that Eq. (@J) 
is written 



R" + k 2 R + eaR 3 + eV Rcos(Kx) = 
and Eq. © is written 
d 2 Ri 



K 2 Ri 



hj(£,R k ,D l R k ) 



(7) 



(8) 



where we recall that n — ex, £ = bx = (1 + eb\ + e 2 b 2 + 
• • • )x, and D l Rk signifies the presence of derivatives of 
i?^ in the right hand side of the equation. Because of the 
scaling in Eq. J7J), ho = 0, so that the 0(1) term is an 
unforced harmonic oscillator. Its solution is 



R = A (r)) cos«) + Bo (77) sin«) . 



(9) 



where Aq and Bq will be determined by the solvability 
condition at O(e). 

The 0(e- ? ) (j > 1) equations arising from Q are forced 
harmonic oscillators, with forcing terms depending on the 
previously obtained (£,??) (k < j) and their deriva- 
tives. Their solutions take the form 



Rj = A 3 {n) cos«) + Bj{rj) sin«) + R 



(10) 



where each Rj p contains contributions from various har- 
monics. As sinusoidal terms giving a 1:1 resonance with 
the OL arise at 0(e 2 ), we can stop at that order. 

At 0(e), there is a contribution from both the OL and 
the nonlinearity, giving 



hx = -V R cos«) - ai?o - 2 



d 2 R 



2b x 



d 2 R 



(11) 



where we recall that the OL depends on the stretched 
spatial variable £ because we are detuning from a reso- 
nant state [23. With Eq. ©, we obtain 



hi = 



2hK 2 A 



2b 1 K 2 B 



2kB' 



aA {Al 



2kA' - -aB (A 2 



oAq 

4 
a Bo 
~ ' 
V0A0 



-3A 



35 2 ) cos(3«£) 
f Bq) sin(3<) - 



B 2 ) 
V A 



cos(k£) 
sin(ft£) 



cos(2k£) 



V B 



sin(2K^) 



(12) 



For rj) to be bounded, the coefficients of the sec- 

ular terms in Eq. 1)12(1 must vanish [2(1 l22l | . The only 
secular harmonics are cos(k£) and sin(K^), and equating 
their coefficients to zero yields the following equations of 
motion describing the slow dynamics: 



A' = -b lK B + ^B (A 2 + B 2 ) , 

cm 

B'^b^Av-^-AoiAl + B 2 ). 



(13) 



We convert l|13fl to polar coordinates with Ao(r]) = 
Cq cos[(/?o( 7 ?)] and Bo(r/) = Co sm[(po(r])] and see imme- 
diately that each circle of constant Co is invariant. The 
dynamics on each circle is given by 

^) = ^(0)+(M-^jr,. (14) 

We examine the special circle of equilibria, corresponding 
to periodic orbits of Q, which satisfies 



Cl = A 2 + B 2 = 



3o: 



(15) 



In choosing an initial configuration for numerical simula- 
tions of the GP equation we set Bo = without loss 
of generality. 

Equating coefficients of © at 0(e 2 ) yields 



d 2 R- 



2 9 

' k R 2 = h 2 



(16) 



where the forcing term again contains contributions from 
both the OL and the nonlinearity: 



d 2 Ro d 2 R 



h 2 = -(bl+2b 2 )^ — ~ -2b, 
d 2 R 1 d 2 R x 



di] 2 



SaRoRi - 2b v 



dS, 2 2 d£dr] 



d 2 R 
d^drj 

— R\Vo cos(k£) . 

(17) 



Here, one inserts the expressions for Ro, Ri, and their 
derivatives into the function h 2 . 

To find the secular terms in Eq. 117(1 . we compute 

Rift, 77) = Ai(r)) cos«) + Si (77) sin«) + R lp (£, 77) , 
Rlp{€> v) = c i cos(3k£) + c 2 sin(3fi;£) 

+ c 3 + C4 cos(2k£) + C5 sin(2re£) , (18) 



where 



Cl 



C3 



—A Q (A 2 ~3B 2 ), C2 = ^B (3A 2 -B 2 ). 



VqAq 
2k 2 



VqAq 
6k 2 



32k 2 
VqBq 
6k 2 



(19) 



After it is expanded, the function h 2 contains harmon- 
ics of the form cos(0£) = 1, cos(k£) (the secular terms), 
cos(2k£), cos(3k£), cos(4k£), and cos(5k£) (as well as 



7 



sine functions with the same arguments). Equating the 
secular cofficients to zeros gives the following equations 
describing the slow dynamics: 



4 



1 



[(/i(a, k)B 2 q + f 2 (a, n)A 2 Q + f 3 (a, re, h)) B 1 



3072k 5 

+/ 4 (a, k)A BqA 1 + f 5 (a, k)B% 
+/ 6 (a, n)A 2 Bl + f 7 (a, k)A^B + / 8s (a, k, b 2 )B ] 
1 



[(/i(a, k)A 2 + fc(a, k)B 2 + f 3 (a, re, 6 X )) A 1 



3072k 5 

+U{a, ^AqBqBx + / 5 (a, re)A;> 



where 



(a, re)^B 2 


+ / 7 (a,K)A S 4 


+ / 8c (o;,re)Ao 


/i(a, re) 


= 3/ 2 (a, «) , 




/2(a, «) 


= -1152are 4 , 




/ 3 (a,/c,&i) 


= 3072k 6 6i , 




/ 4 (a, re) 


= 2/2 (a, re) , 




/s(a, re) 


= 180a 2 re 2 , 




/e(a, re) 


= 2/5 (a, re) , 




/ 7 (a, re) 


= /5(a,re) , 




f 8s (a,K,b 2 ) 


= f non (a, re) - 


128V r 2 re 2 , 


/8c (a, re) 


= fnon(a, re) + 


640 V 2 re 2 , 


fnon if^z ^) 


= 3072k 6 & 2 ■ 





(20) 



(21) 



We use the notation f non to indicate the portions of the 
quantities fg s and / 8c that arise from non-resonant terms. 
The other terms in these quantities, which depend on the 
lattice amplitude Vq, arise from resonant terms. 



Equilibrium solutions of l|20|) satisfy 



(fiBj + hA 2 + h){hAl + UAlB 2 + f 7 A Bl + / 8c A ) - (f 4 A B )(f 5 Bl + / 6 A 2 gg + / 7 A 4 5 + f 8s B ) 

flA 2 Bl - {hBl + f 2 A 2 + h){hAl + frBl +h) 
{hA 2 + / 2 g 2 + / 3 )(/ 5 gg + hA 2 Bl + f 7 A*B + f Bs B ) - (UA B )(f 5 A 5 Q + f 6 A 3 B 2 + f 7 A B 4 + / 8e A ) 

/| A\Bl {hBl + f 2 A 2 + h){hA 2 + hBl +h) 

(22) 



where one uses an equilibrium value of Aq and Bq from 
Eq. i|15|) . Inserting equilibrium values of Aq, Bq, Ai, and 
i?! into Eqs. 10 and (|18fl. we obtain the spatial profile 
R = Rq + eRi + 0(e 2 ) used as the initial wavefunction 
in the numerical simulations of the full GP equation JTJ 
with a stationary OL. 



VI. APPENDIX B: ANALYTICAL 
CONSTRUCTION OF QUARTER-PERIOD 
STATES 

To construct quarter-period states, we use the reso- 
nance relation yj~8 = 2k and the scaling a = £ 4 a, so that 



I 

Eq. J2J is written 

R" + Ak 2 R + s 4 aR 3 + sVqR cos(kx) = (23) 
and Eq. © is written 

Lt[Rj] = ^ + 4k 2 E, = fitf, Rk,D l R k ) , (24) 

where 77 = ex and £ = bx = (1 + eb\ + e 2 b 2 + ■ ■ ■ )x, as 
before. 

Because of the scaling in (I23f) . ho = (as in the case of 
half-period states), so that the 0(1) term is an unforced 
harmonic oscillator. It has the solution 

R = A {rj) cos(2kcI) + B { v ) sin(2re0 , (25) 
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where Aq 



B 2 



Cq is an arbitrary constant (in the where 



numerical simulations, we take B — without loss of 
generality) . With the different scaling of the nonlinearity 
coefficient, the value C 2 is not constrained as it was in 
the case of half-period states (see Appendix A). 

The 0(e J ) (J > 1) equations arising from (|24|l are 
forced harmonic oscillators, with forcing terms depend- 
ing on the previously obtained i?fc(£, Tj) (k < j) and their 
derivatives. Their solutions take the form 



Rj = Aj(rj) cos(2k£) + Bj(r)) sui(2k£) 



(26) 



where Rj P contain contributions from various harmonics. 
As sinusoidal terms giving a 1:2 resonance with the OL 
arise at 0(e ), we can stop at that order. 

The equation at 0(e) has a solution of the form 



Ri = Ai(r)) cos(2<) + B^rf) sm(2fi£) + R 



(27) 



The coefficients A\ and B\ are determined using a solv- 
ability condition obtained at 0(e 2 ) by requiring that the 
secular terms of h 2 vanish. This yields 



At 

Bi = 



Vn 



I661K 2 
% 



I661K 2 

The particular solution is 



(en 
(C13 



C12) - — A 
bi 



c M ) 



bi 



-Bo 



(28) 



-Rip = en cos(k£ )+c\2 cos(3k^)+ci3 sin(3K^)+Ci4 sui(k£) , 

(29) 

where 



Cll 



Cl3 



V Ao 


V A 


6k 2 


10n z 


V B 


VqBq 


10k 2 ' 


Cl4= 6k 2 



(30) 



The solution at 0{e 2 ) has the form 

R2 = A 2 (r])cos(2K^+B 2 (r])sm(2KO + R2p^,v)- (31) 

The coefficients A 2 and B 2 are determined using a solv- 
ability condition obtained at 0(e 3 ) by requiring that the 
secular terms of h 3 vanish. This yields 



A* 



B> 



Vo 



I661K 2 
b 3 



(C21 + C22) 



An 



Vo 
32k 2 

b_2 V 2 A 
bi 1 480k 4 : 



(en -I- C12) 



^0 



16felK 2 



(C23 + C24) 



Vo 



(Cl 3 + C14) 



63 bo 
- 7^0 - -^B 1 



V 2 B 



480k 4 ■ 
The particular solution is 

R2 P = C21 cos(k£) + c 22 cos(3k£) 
+ c 2 4 sui(k£) + C25 cos(4k£) 
+ c 2 6 + C27 sin(4K^) , 



(32) 



c 23 sui(3k£) 



(33) 







h^ Vn /In 


C21 = 


6k 2 


9K 




VqAi 




C22 = 


10k 2 


25k 2 ' 




VoBi 


3&lFo-Bo 


C23 = 


10k 2 


25k 2 




%Bt 


hV B 


C24 = 


6k 2 


9k 2 



C25 



240k 4 



C26 



48k 4 



C27 



Vo 2 Bq 
240k 4 



(34) 



Note that the harmonics cos(0£) and sin(0£) occur in 133|) 
and are reduced appropriately. (The arguments of this 
sine and cosine arise because of our particular resonance 
relation.) 

At 0(e 3 ), we obtain solutions of the form 

R 3 = A 3 ( v ) cos(2<) + B 3 ( v ) sin(2<) + R 3p (^, v) • (35) 

The coefficients A 3 and B 3 are determined using a solv- 
ability condition obtained at 0(e 4 ) by requiring that the 
secular terms of /14 vanish. Because of the scaling in l|23|l . 
the effects of the nonlinearity manifest in this solvability 
condition. The resulting coefficients are 



A 3 = - 



Vq% b 3 V 2 b 2 A 
33756ik 8 61 1 1800&1K 4 
V 2 Ai miVgAo daAf, 



1800k 4 
b2 A 



54000k 4 326ik 2 

V 2 A 2 , iaA B 2 



bi 



A n 



2406ik 4 326ik 2 



(36) 



B 3 = 



b 2 

-r~B 2 -+ 

1800k 4 61 1 
V 2 b 2 B 



3aAlB Q 
32b lK 2 



64 B 
7- Bo 

b\ 



VqB 2 
2406ik 4 



119V^B 



8640006ik 8 
l9hV 2 B 



ZaBl 



I8OO61K 4 
The particular solution is 

R^p 



54000k 4 32&ik 2 



(37) 



c 3 i cos(k£) + c 32 cos(3k£) + C33 sin(3K£) 
C34 sm(K^) + C35 cos(4k£) + c 36 
C37 sin(4K£) -I- c 38 cos(5k£) + c 39 cos(k£) 
C310 sui(5k£) + C311 sin(K£) , (38) 



where 



C31 = 



7b 2 V A , b 2 V Ao 11V 3 A 



54k 2 
9k 2 



9k 2 
VqA 2 
6k 2 



4320k 6 



(39) 



9 



C32 



C33 



mlVpAo Sb 2 V Q A Q 17V 3 A 
250k 2 25k 2 12000k 6 



25k 2 



25k 2 
VqA 2 
10k 2 ' 



17V S B 3&iV Si ^ 316fV B 
12000k 6 25k 2 ^ 



25k 2 



10k 2 



250k 2 



Similar to what occurs at 0(e 2 ), the coefficient C36 is 
the prefactor for cos(0£) and a sin(0£) term (not shown) 
occurs in (|38|l as well. The extra terms (from the reso- 

(40) 

nance relation) that go into the slow evolution equations 
and the resulting expressions for the periodic orbits (i.e., 
the equilibria of the slow flow) arise from the terms with 
prefactors C39 and C311. (The harmonics corresponding 
to the coefficients C31 and C34 are always secular, but 
those corresponding to C39 and C311 are secular only for 
(41) 1:2 superharmonic resonances.) 



f'34 



hVoBp lb\V Q B V B 2 

9k 2 
hVoBt 



9k 2 



54k 2 
4320k 6 



6k 2 



(42) 



C35 
C37 
C39 



19&iV 2 4> V 2 A 1 



1800k 4 
19feiVo 2 go 

1800k 4 
V 3 A 



288k 6 



C310 



240k 4 

240k 4 

% 3 B 



C36 



C-38 



via l _ 

48k 4 

VjAo 
10080k 6 



biV 2 A 
72k 4 



10080k 6 



C311 = 



Vq 3 B 
288k 6 



(43) 



One inserts equilibrium values of Aj and Bj (j € 
{0,1,2,3}) into Eqs. (gSJl, 0, JSD), and ESJ to obtain 
the spatial profile i? = i? + e-Ri + £ 2 B-2 +e 3 R 3 + 0(e 4 ) 
used as the initial wavefunction in numerical simulations 
of the GP equation with a stationary OL. 



[1] S. Aubry, Physica D 103, 201, (1997); S. Flach and C. 
R. Willis, Phys. Rep. 295, 181 (1998); D. Hennig and G. 
Tsironis, Phys. Rep. 307, 333 (1999); P. G. Kevrekidis, 
K. O. Rasmussen, and A. R. Bishop, Int. J. Mod. Phys. 
B 15, 2833 (2001). 

[2] D. K. Campbell, P. Rosenau, and G. Zaslavsky, Chaos 
15, 015101 (2005). 

[3] D. N. Christodoulides, F. Lederer, and Y. Silberberg, Na- 
ture 424, 817 (2003); Yu. S. Kivshar and G. P. Agrawal, 
Optical Solitons: From Fibers to Photonic Crystals (Aca- 
demic Press, San Diego, 2003); J. W. Fleischer, G. Bar- 
tal, O. Cohen, T. Schwartz, O. Manela, B. Freedman, M. 
Segev, H. Buljan, and N. K. Efremidis, Opt. Express 13, 
1780 (2005). 

[4] P. G. Kevrekidis and D. J. Frantzeskakis, Mod. Phys. 
Lett. B 18, 173 (2004); V. V. Konotop and V. A. Brazh- 
nyi, Mod. Phys. Lett. B 18 627, (2004); P.G Kevrekidis, 
R. Carretero-Gonzalez, D. J. Frantzeskakis, and I. G. 
Kevrekidis, Mod. Phys. Lett. B 18, 1481 (2004); O. 
Morsch and M. Oberthaler, Rev. Mod. Phys. 78, 179 
(2006). 

[5] M. Peyrard, Nonlinearity 17, Rl (2004). 

[6] A. Trombettoni and A. Smerzi, Phys. Rev. Lett. 86, 2353 

(2001) ; F. Kh. Abdullaev, B. B. Baizakov, S. A. Dar- 
manyan, V. V. Konotop, and M. Salerno, Phys. Rev. A 
64, 043606 (2001); G. L. Alfimov, P. G. Kevrekidis, V. 
V. Konotop, and M. Salerno, Phys. Rev. E 66, 046608 

(2002) . 

[7] J. P. Keener, Phys. D 136, 1 (2000); P. G. Kevrekidis 

and I. G. Kevrekidis, Phys. Rev. E 64, 056624 (2001). 
[8] T. Stoferle, H. Moritz, C. Schori, M. Kohl, and T. 



Esslinger, Phys. Rev. Lett. 92, 130403 (2004). 
[9] M. Kramer, C. Tozzo, and F. Dalfovo, Phys. Rev. A 71, 
061602(R) (2005). 
[10] C. Tozzo, M. Kramer, and F. Dalfovo, Phys. Rev. A 72, 
023613 (2005). 

[11] A. Iucci, M. A. Cazalilla, A. F. Ho, and T. Giamarchi, 

Phys. Rev. A. 73, 041608 (2006). 
[12] N. Gemelke, E. Sarajlic, Y. Bidel, S. Hong, and S. Chu, 

Phys. Rev. Lett. 95, 170404 (2005). 
[13] M. Machholm, A. Nicolin, C. J. Pethick, and H. Smith, 

Phys. Rev. A 69, 043604 (2004). 
[14] M. A. Porter and P. Cvitanovic, Phys. Rev. E 69, 047201 

(2004); M. A. Porter and P. Cvitanovic, Chaos 14, 739 

(2004). 

[15] A. Smerzi, A. Trombettoni, P. G. Kevrekidis, and A. R. 
Bishop, Phys. Rev. Lett. 89, 170402 (2002). 

[16] F. S. Cataliotti, L. Fallani, F. Ferlaino, C. Fort, P. Mad- 
daloni, and M. Inguscio, New J. Phys. 5, 71 (2003). 

[17] B. T. Seaman, L. D. Carr, and M. J. Holland, Phys. Rev. 
A 72, 033602 (2005). 

[18] N. W. Ashcroft and N. D. Mermin, Solid State Physics 
(Saunders College, Philadelphia, 1976). 

[19] V. M. Perez-Garcfa, H. Michinel, and H. Herrero, Phys. 
Rev. A 57, 3837 (1998); L. Salasnich, A. Parola, and L. 
Reatto, Phys. Rev. A 65, 043614 (2002). 

[20] A. H. Nayfeh and D. T. Mook, Nonlinear Oscillations 
(John Wiley & Sons, New York, 1995). 

[21] The small sizes of the lattice depth and nonlinearity co- 
efficient allow us to perturb from base states consisting 
of trigonometric functions (regular harmonics). One can 
relax these conditions using the much more technically 



10 



demanding approaches of perturbing from either Math- 
ieu functions (removing the restriction on the size of the 
amplitude of the periodic potential) or elliptic functions 
(removing the restriction on the size of the nonlinearity 
coefficient). 

[22] R. H. Rand, Lecture Notes on Nonlin- 
ear V ibrations, online book available at 
\http://www. tam. Cornell. edu/randdocs/nlvibe52.pdj\ 
2005. 

[23] M. A. Porter and P. G. Kevrekidis, SIAM J. App. Dyn. 

Sys. 4, 783 (2005). 
[24] N. W. McLachlan, Theory and Application of Mathieu 



Functions (Clarendon Press, Oxford, 1947). 

[25] We use the same values of Aq and Bo (the prefactors 
of the harmonics in _Ro; see the appendices) as for the 
half-period states, even though these parameters are not 
constrained for quarter-period solutions. As discussed in 
Appendix A, they are constrained for half-period states. 

[26] In many cases, one can make this procedure more math- 
ematically rigorous (though less transparent physically) 
by examining the dynamics geometrically and introduc- 
ing slow and fast manifolds. 



